microphone differences?

## subject & microphone list ----

## fomri arrive (range):
fomri.arrive.lb <- "2018-05-01"
fomri.arrive.ub <- "2018-05-31"

## definite fomri
fomri.switch <- "2018-06-28"
## exceptions:
##  DMCC515268_baseline DMCC3 (June 30th, 2018) and DMCC6904377_reactive DMCC2 (July 1st, 2018)

dates1.is.microoptics <- dates1[c("BAS.day", "PRO.day", "REA.day")] < fomri.arrive.lb
dates1.is.fomri <- dates1[c("BAS.day", "PRO.day", "REA.day")] > fomri.switch

all(dates1.is.microoptics, na.rm = TRUE)  ## all microoptics
## [1] TRUE
dates2.is.microoptics <- dates2[c("session1", "session2", "session3")] < fomri.arrive.lb
dates2.is.fomri <- dates2[c("session1", "session2", "session3")] > fomri.switch

dates2$all.microoptics <- apply(dates2.is.microoptics, 1, function(x) all(x))
dates2$all.fomri <- apply(dates2.is.fomri, 1, function(x) all(x))

subj.microoptics <- dates2 %>% filter(all.microoptics & subj != "DMCC6904377") %>% .$subj
subj.microoptics <- union(subj.microoptics, dates1$subj)
subj.fomri <- dates2 %>% filter(all.fomri & subj != "DMCC6904377") %>% .$subj


stroop$mic <- ifelse(
  stroop$subj %in% subj.fomri,
  "fomri",
  ifelse(
    stroop$subj %in% subj.microoptics, 
    "micro",
    "unknown"
  )
)

table(stroop$subj, stroop$mic)
##              
##               fomri micro unknown
##   107321          0     0     672
##   115825          0   672       0
##   123117          0   672       0
##   130114          0   672       0
##   130518          0   672       0
##   132017          0   672       0
##   138837          0   672       0
##   141422          0   672       0
##   150423          0   672       0
##   158136          0   672       0
##   160830          0     0     672
##   161832          0     0     672
##   165032          0   672       0
##   171330          0   456       0
##   178243          0   645       0
##   178647          0   672       0
##   178950          0   672       0
##   182436          0     0     672
##   197449          0   672       0
##   203418          0   672       0
##   204319          0   672       0
##   233326          0   672       0
##   300618          0   672       0
##   317332          0   672       0
##   346945          0   672       0
##   352738          0   432       0
##   393550          0   672       0
##   448347        672     0       0
##   580650          0   672       0
##   594156          0   672       0
##   601127          0   672       0
##   765864          0   672       0
##   814649          0   672       0
##   843151          0   456       0
##   849971          0   672       0
##   873968        672     0       0
##   877168          0   672       0
##   DMCC1328342   672     0       0
##   DMCC1596165   672     0       0
##   DMCC1971064     0   672       0
##   DMCC2442951   672     0       0
##   DMCC2803654     0   672       0
##   DMCC2834766     0   672       0
##   DMCC3062542   672     0       0
##   DMCC4191255     0   672       0
##   DMCC4260551     0     0     672
##   DMCC4854075     0   672       0
##   DMCC5009144   456     0       0
##   DMCC5195268     0   672       0
##   DMCC5775387     0   672       0
##   DMCC6371570   672     0       0
##   DMCC6418065     0   672       0
##   DMCC6484785     0     0     672
##   DMCC6661074     0   432       0
##   DMCC6671683     0   672       0
##   DMCC6705371     0   672       0
##   DMCC6721369     0   672       0
##   DMCC6904377     0     0     456
##   DMCC7297690   432     0       0
##   DMCC7921988     0     0     672
##   DMCC8033964     0   456       0
##   DMCC8050964   672     0       0
##   DMCC8078683     0     0     672
##   DMCC9441378   672     0       0
##   DMCC9478705     0   672       0
##   DMCC9953810   672     0       0
stroop %>%
  ggplot(aes(subj, rt, color = mic)) +
  facet_grid(vars(session)) +
  geom_point(position = position_jitter(width = 0.1), size = 0.5, alpha = 0.3) +
  geom_boxplot(position = position_nudge(1/3), notch = TRUE, width = 0.2) +
  scale_color_brewer(type = "qual", palette = 2) +
  theme(axis.text.x = element_text(angle = 90, hjust = 0, vjust = 0))
## Warning: Removed 422 rows containing non-finite values (stat_boxplot).
## Warning: Removed 422 rows containing missing values (geom_point).
stroop %>%
  group_by(mic, subj) %>%
  summarize(freq = sum(rt == 0, na.rm = TRUE)) %>%
  ggplot(aes(mic, freq)) +
  geom_point(size = 2) +
  labs(y = "frequency of rt == 0 per subj")
stroop %>%
  group_by(mic, subj, error) %>%
  summarize(freq = sum(error, na.rm = TRUE)) %>%
  ggplot(aes(mic, freq)) +
  geom_point(size = 2, position = position_jitter(width = 0.1), alpha = 0.4) +
  labs(y = "frequency of errors per subj")
stroop %>%
  group_by(mic, subj, acc.final) %>%
  summarize(freq = n()) %>%
  ggplot(aes(mic, freq)) +
  facet_grid(cols = vars(acc.final)) +
  geom_point(size = 2, position = position_jitter(width = 0.1), alpha = 0.4) +
  labs(y = "frequency of acc.final codings per subj")
stroop %>%
  ggplot(aes(sample = rt, color = mic)) +
  stat_qq(alpha = 0.3, size = 0.15) +
  stat_qq_line(size = 0.25) +
  facet_wrap(vars(subj)) +
  scale_color_brewer(type = "qual", palette = 2)
## Warning: Removed 422 rows containing non-finite values (stat_qq).
## Warning: Removed 422 rows containing non-finite values (stat_qq_line).
stroop.mean.sd <- stroop %>%
  group_by(mic, subj, session) %>%
  summarize(rt.mean = mean(rt, na.rm = TRUE), rt.sd = sd(rt, na.rm = TRUE))

grid.arrange(
  stroop.mean.sd %>%
    ggplot(aes(mic, rt.mean)) +
    facet_grid(vars(session)) +
    geom_point(position = position_jitter(width = 0.1)) +
    geom_boxplot(position = position_nudge(1/3), notch = FALSE, width = 0.2),
  stroop.mean.sd %>%
    ggplot(aes(mic, rt.sd)) +
    facet_grid(vars(session)) +
    geom_point(position = position_jitter(width = 0.1)) +
    geom_boxplot(position = position_nudge(1/3), notch = FALSE, width = 0.2)
)

model

differences in mean stroop effect btw mics?

stroop.rt <- stroop %>% filter(!is.na(rt), rt > 200, acc == 1)

# g <- paste0(stroop.rt$subj, "_", stroop.rt$session)
# l <- split(stroop.rt, g)
# 
# fits <- lapply(l, function(d) lm(rt ~ trial.type, d))
# stroop.rt$resids <- unlist(lapply(fits, resid), use.names = FALSE)
# bs <- do.call(rbind, lapply(fits, coef))
# bs <- cbind(bs, reshape2::colsplit(rownames(bs), "_", c("subj", "session"))

stroop.rt %<>%
  mutate(
    bas   = as.numeric(session == "bas"),
    pro   = as.numeric(session == "pro"),
    rea   = as.numeric(session == "rea"),
    congr = as.numeric(trial.type == "c"),
    incon = as.numeric(trial.type == "i"),
    micro = as.numeric(mic == "micro"),
    fomri = as.numeric(mic == "fomri")
  )

fit.noint <- lmer(
  rt ~ 
    0 +
    congr:bas:micro + congr:pro:micro + congr:rea:micro +
    incon:bas:micro + incon:pro:micro + incon:rea:micro +
    congr:bas:fomri + congr:pro:fomri + congr:rea:fomri +
    incon:bas:fomri + incon:pro:fomri + incon:rea:fomri +
    (
    0 + 
    congr:bas + congr:pro + congr:rea +
    incon:bas + incon:pro + incon:rea
    | subj
    ),
  stroop.rt,
  subset = mic != "unknown",
  control = lmerControl(optimizer = "bobyqa")
  )

# fits.noint <- allFit(fit.noint)
# summary(fits.noint)
# fit.noint <- fits.noint$bobyqa

summary(fit.noint)
## Linear mixed model fit by REML ['lmerMod']
## Formula: rt ~ 0 + congr:bas:micro + congr:pro:micro + congr:rea:micro +  
##     incon:bas:micro + incon:pro:micro + incon:rea:micro + congr:bas:fomri +  
##     congr:pro:fomri + congr:rea:fomri + incon:bas:fomri + incon:pro:fomri +  
##     incon:rea:fomri + (0 + congr:bas + congr:pro + congr:rea +  
##     incon:bas + incon:pro + incon:rea | subj)
##    Data: stroop.rt
## Control: lmerControl(optimizer = "bobyqa")
##  Subset: mic != "unknown"
## 
## REML criterion at convergence: 467003
## 
## Scaled residuals: 
##    Min     1Q Median     3Q    Max 
## -5.019 -0.532 -0.123  0.347 12.417 
## 
## Random effects:
##  Groups   Name      Variance Std.Dev. Corr                    
##  subj     congr:bas 14727    121                              
##           congr:pro 12694    113      0.71                    
##           congr:rea 18702    137      0.71 0.84               
##           bas:incon 19455    139      0.89 0.66 0.60          
##           pro:incon 16440    128      0.72 0.97 0.79 0.76     
##           rea:incon 20603    144      0.75 0.87 0.94 0.75 0.90
##  Residual           27957    167                              
## Number of obs: 35655, groups:  subj, 57
## 
## Fixed effects:
##                 Estimate Std. Error t value
## congr:bas:micro    830.0       18.5    44.9
## congr:micro:pro    821.6       17.1    48.2
## congr:micro:rea    810.7       20.6    39.3
## bas:micro:incon    962.4       21.3    45.2
## micro:pro:incon    900.2       19.2    46.8
## micro:rea:incon    909.4       21.6    42.0
## congr:bas:fomri    734.6       36.0    20.4
## congr:pro:fomri    718.8       33.0    21.8
## congr:rea:fomri    767.1       40.2    19.1
## bas:incon:fomri    860.2       41.3    20.8
## pro:incon:fomri    810.0       37.2    21.8
## rea:incon:fomri    847.6       42.1    20.1
## 
## Correlation of Fixed Effects:
##             cngr:bs:m cngr:mcr:p cngr:mcr:r bs:mc: mcr:p: mcr:r: cngr:bs:f
## cngr:mcr:pr 0.684                                                         
## congr:mcr:r 0.689     0.813                                               
## bas:mcr:ncn 0.878     0.635      0.579                                    
## micr:pr:ncn 0.698     0.945      0.774      0.737                         
## micro:r:ncn 0.724     0.851      0.923      0.728  0.883                  
## cngr:bs:fmr 0.000     0.000      0.000      0.000  0.000  0.000           
## cngr:pr:fmr 0.000     0.000      0.000      0.000  0.000  0.000  0.682    
## congr:r:fmr 0.000     0.000      0.000      0.000  0.000  0.000  0.681    
## bas:ncn:fmr 0.000     0.000      0.000      0.000  0.000  0.000  0.878    
## pro:ncn:fmr 0.000     0.000      0.000      0.000  0.000  0.000  0.696    
## rea:ncn:fmr 0.000     0.000      0.000      0.000  0.000  0.000  0.718    
##             cngr:p: cngr:r: bs:nc: pr:nc:
## cngr:mcr:pr                              
## congr:mcr:r                              
## bas:mcr:ncn                              
## micr:pr:ncn                              
## micro:r:ncn                              
## cngr:bs:fmr                              
## cngr:pr:fmr                              
## congr:r:fmr 0.808                        
## bas:ncn:fmr 0.634   0.574                
## pro:ncn:fmr 0.946   0.769   0.734        
## rea:ncn:fmr 0.848   0.922   0.723  0.879
bnames <- rownames(coef(summary(fit.noint)))

m <- data.frame(
  bas   = grepl("bas", bnames),
  pro   = grepl("pro", bnames),
  rea   = grepl("rea", bnames),
  congr = grepl("congr", bnames),
  incon = grepl("incon", bnames),
  fomri = grepl("fomri", bnames),
  micro = grepl("micro", bnames)
)

rownames(m) <- bnames
contrast <- rbind(
  "fomri-micro|bas" = (m$fomri - m$micro) * m$bas / 2,
  "fomri-micro|pro" = (m$fomri - m$micro) * m$pro / 2,
  "fomri-micro|rea" = (m$fomri - m$micro) * m$rea / 2,
  "stroop:(fomri-micro)|bas" = (m$incon - m$congr) * (m$fomri - m$micro) * m$bas / 2,
  "stroop:(fomri-micro)|pro" = (m$incon - m$congr) * (m$fomri - m$micro) * m$pro / 2,
  "stroop:(fomri-micro)|rea" = (m$incon - m$congr) * (m$fomri - m$micro) * m$rea / 2
)
colnames(contrast) <- bnames

results <- glht(fit.noint, contrast, alternative = "two.sided")
summary(results, adjusted(type = "none"))
## 
##   Simultaneous Tests for General Linear Hypotheses
## 
## Fit: lmer(formula = rt ~ 0 + congr:bas:micro + congr:pro:micro + congr:rea:micro + 
##     incon:bas:micro + incon:pro:micro + incon:rea:micro + congr:bas:fomri + 
##     congr:pro:fomri + congr:rea:fomri + incon:bas:fomri + incon:pro:fomri + 
##     incon:rea:fomri + (0 + congr:bas + congr:pro + congr:rea + 
##     incon:bas + incon:pro + incon:rea | subj), data = stroop.rt, 
##     control = lmerControl(optimizer = "bobyqa"), subset = mic != 
##         "unknown")
## 
## Linear Hypotheses:
##                               Estimate Std. Error z value Pr(>|z|)  
## fomri-micro|bas == 0            -98.83      42.12   -2.35    0.019 *
## fomri-micro|pro == 0            -96.56      39.00   -2.48    0.013 *
## fomri-micro|rea == 0            -52.71      45.35   -1.16    0.245  
## stroop:(fomri-micro)|bas == 0    -3.36      11.12   -0.30    0.762  
## stroop:(fomri-micro)|pro == 0     6.27       6.93    0.91    0.365  
## stroop:(fomri-micro)|rea == 0    -9.08       9.18   -0.99    0.322  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## (Adjusted p values reported -- none method)

Summary

  • fomri microphone recorded faster RTs
  • no observed impact of microphone on size of mean stroop effect

differences in variance btw microphones?

cl1 <- lmeControl(
 maxIter = 100000, msMaxIter = 100000, niterEM = 100000,
 msMaxEval = 100000, tolerance = 0.000001, msTol = 0.0000001, returnObject = TRUE,
 minAbsParApVar = 0.05, opt = c("nlminb"), optimMethod = "BFGS"
 )
cl2 <- lmeControl(
 maxIter = 100000, msMaxIter = 100000, niterEM = 100000,
 msMaxEval = 100000, tolerance = 0.000001, msTol = 0.0000001, returnObject = TRUE,
 minAbsParApVar = 0.05, opt = c("optim"), optimMethod = "BFGS"
 )
cl3 <- lmeControl(
 maxIter = 100000, msMaxIter = 100000, niterEM = 100000,
 msMaxEval = 100000, tolerance = 0.000001, msTol = 0.0000001, returnObject = TRUE,
 minAbsParApVar = 0.05, opt = c("optim"), optimMethod = "L-BFGS-B"
 )

## proactive

m.hom.pro <- lme(
  rt ~ trial.type * mic, 
  random = ~ 1 | subj,
  data = stroop.rt %>% filter(mic != "unknown"), 
  subset = session == "pro",
  control = cl1
  )

m.het.mic.pro <- update(m.hom.pro, . ~ ., weights = varIdent(form = ∼ 1 | mic))
m.het.subj.pro <- update(m.hom.pro, . ~ ., weights = varIdent(form = ∼ 1 | subj))

summary(m.hom.pro)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "pro" 
##      AIC    BIC logLik
##   157167 157211 -78578
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:         122    169.6
## 
## Fixed effects: rt ~ trial.type * mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          718.8     35.69 11915  20.140  0.0000
## trial.typei           91.1      7.09 11915  12.839  0.0000
## micmicro             102.7     40.17    55   2.557  0.0134
## trial.typei:micmicro -12.4      8.00 11915  -1.546  0.1220
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.132              
## micmicro             -0.888  0.117       
## trial.typei:micmicro  0.117 -0.887 -0.132
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -4.1408 -0.5416 -0.1393  0.3415 12.1557 
## 
## Number of Observations: 11974
## Number of Groups: 57
summary(m.het.mic.pro)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "pro" 
##      AIC    BIC logLik
##   157121 157173 -78553
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:         122    173.4
## 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | mic 
##  Parameter estimates:
##  micro  fomri 
## 1.0000 0.8939 
## Fixed effects: rt ~ trial.type + mic + trial.type:mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          718.8     35.61 11915  20.183  0.0000
## trial.typei           91.1      6.48 11915  14.052  0.0000
## micmicro             102.7     40.11    55   2.561  0.0132
## trial.typei:micmicro -12.4      7.50 11915  -1.649  0.0992
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.121              
## micmicro             -0.888  0.107       
## trial.typei:micmicro  0.104 -0.864 -0.124
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -4.0508 -0.5440 -0.1399  0.3416 11.8919 
## 
## Number of Observations: 11974
## Number of Groups: 57
summary(m.het.subj.pro)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "pro" 
##      AIC    BIC logLik
##   153344 153802 -76610
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:       121.8    159.9
## 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | subj 
##  Parameter estimates:
##      115825      123117      130114      130518      132017      138837 
##      1.0000      1.3684      0.9483      0.7232      1.1539      1.0273 
##      141422      150423      158136      165032      171330      178243 
##      1.0303      1.3350      0.6193      1.6752      0.9424      1.5864 
##      178647      178950      197449      203418      204319      233326 
##      1.1050      0.7494      1.1987      0.7306      1.8857      2.0556 
##      300618      317332      346945      352738      393550      448347 
##      0.5359      0.6395      0.8655      1.6850      0.7852      1.1969 
##      580650      594156      601127      765864      814649      843151 
##      2.0928      0.6314      1.1921      0.5608      0.9475      1.1984 
##      849971      873968      877168 DMCC1328342 DMCC1596165 DMCC1971064 
##      1.8200      0.7821      1.3117      0.4071      1.6184      0.8052 
## DMCC2442951 DMCC2803654 DMCC2834766 DMCC3062542 DMCC4191255 DMCC4854075 
##      0.4382      0.6706      0.7042      0.9748      0.6480      0.6854 
## DMCC5009144 DMCC5195268 DMCC5775387 DMCC6371570 DMCC6418065 DMCC6661074 
##      0.5089      0.5855      1.3942      0.9953      0.5380      0.7065 
## DMCC6671683 DMCC6705371 DMCC6721369 DMCC7297690 DMCC8033964 DMCC8050964 
##      0.5219      0.7398      0.5582      0.7739      0.5176      0.8544 
## DMCC9441378 DMCC9478705 DMCC9953810 
##      1.4022      0.8540      0.9206 
## Fixed effects: rt ~ trial.type + mic + trial.type:mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          735.5     35.44 11915  20.756  0.0000
## trial.typei           65.2      4.67 11915  13.940  0.0000
## micmicro              88.6     39.90    55   2.221  0.0305
## trial.typei:micmicro   8.7      5.45 11915   1.589  0.1121
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.087              
## micmicro             -0.888  0.078       
## trial.typei:micmicro  0.075 -0.858 -0.091
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -4.9122 -0.6275 -0.1734  0.4266 10.1056 
## 
## Number of Observations: 11974
## Number of Groups: 57
anova(m.hom.pro, m.het.mic.pro, m.het.subj.pro)
##                Model df    AIC    BIC logLik   Test L.Ratio p-value
## m.hom.pro          1  6 157167 157211 -78578                       
## m.het.mic.pro      2  7 157121 157173 -78553 1 vs 2      48  <.0001
## m.het.subj.pro     3 62 153344 153802 -76610 2 vs 3    3887  <.0001
## baseline

m.hom.bas <- lme(
  rt ~ trial.type * mic, 
  random = ~ 1 | subj,
  data = stroop.rt %>% filter(mic != "unknown"), 
  subset = session == "bas",
  control = cl1
  )

m.het.mic.bas <- update(m.hom.bas, . ~ ., weights = varIdent(form = ∼ 1 | mic))
m.het.subj.bas <- update(m.hom.bas, . ~ ., weights = varIdent(form = ∼ 1 | subj))

summary(m.hom.bas)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "bas" 
##      AIC    BIC logLik
##   144451 144495 -72220
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:       123.1    155.2
## 
## Fixed effects: rt ~ trial.type * mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          742.0     37.31 11101  19.889  0.0000
## trial.typei          130.4      6.80 11101  19.186  0.0000
## micmicro              86.1     41.91    51   2.054  0.0451
## trial.typei:micmicro   3.1      7.66 11101   0.409  0.6829
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.060              
## micmicro             -0.890  0.053       
## trial.typei:micmicro  0.053 -0.888 -0.060
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -4.9472 -0.5558 -0.1074  0.4048 11.8554 
## 
## Number of Observations: 11156
## Number of Groups: 53
summary(m.het.mic.bas)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "bas" 
##      AIC    BIC logLik
##   144351 144402 -72168
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:       123.1    160.2
## 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | mic 
##  Parameter estimates:
##  micro  fomri 
## 1.0000 0.8417 
## Fixed effects: rt ~ trial.type + mic + trial.type:mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          742.0     37.26 11101  19.917  0.0000
## trial.typei          130.4      5.91 11101  22.083  0.0000
## micmicro              86.1     41.87    51   2.056  0.0449
## trial.typei:micmicro   3.1      6.94 11101   0.451  0.6521
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.052              
## micmicro             -0.890  0.047       
## trial.typei:micmicro  0.044 -0.851 -0.054
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -4.7921 -0.5606 -0.1103  0.4044 11.4845 
## 
## Number of Observations: 11156
## Number of Groups: 53
summary(m.het.subj.bas)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "bas" 
##      AIC    BIC logLik
##   142054 142478 -70969
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:       122.4    180.7
## 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | subj 
##  Parameter estimates:
##      115825      123117      130114      130518      132017      138837 
##      1.0000      1.1149      0.9185      0.5873      0.8056      1.1128 
##      141422      150423      158136      165032      178243      178647 
##      0.7795      0.8202      0.6938      0.6761      1.5863      1.1238 
##      178950      197449      203418      204319      233326      300618 
##      0.7375      1.7405      0.8050      0.8296      1.5524      0.6524 
##      317332      346945      352738      393550      448347      580650 
##      0.5346      1.0779      0.6380      0.7650      1.1714      0.6154 
##      594156      601127      765864      814649      849971      873968 
##      0.5954      1.1536      0.5274      1.0802      1.3437      0.7165 
##      877168 DMCC1328342 DMCC1596165 DMCC1971064 DMCC2442951 DMCC2803654 
##      1.1406      0.4654      0.6111      0.8529      0.3401      0.7688 
## DMCC2834766 DMCC3062542 DMCC4191255 DMCC4854075 DMCC5195268 DMCC5775387 
##      0.6989      0.6309      0.5408      0.3850      0.7070      0.7052 
## DMCC6371570 DMCC6418065 DMCC6661074 DMCC6671683 DMCC6705371 DMCC6721369 
##      0.7652      0.7075      0.7614      0.5515      0.5280      0.6040 
## DMCC7297690 DMCC8050964 DMCC9441378 DMCC9478705 DMCC9953810 
##      0.6113      0.9713      0.9493      0.8234      0.6086 
## Fixed effects: rt ~ trial.type + mic + trial.type:mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          747.9     37.03 11101  20.196  0.0000
## trial.typei          112.3      4.76 11101  23.604  0.0000
## micmicro              82.6     41.61    51   1.985  0.0526
## trial.typei:micmicro  12.2      5.60 11101   2.175  0.0297
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.042              
## micmicro             -0.890  0.038       
## trial.typei:micmicro  0.036 -0.849 -0.044
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -4.4068 -0.6214 -0.1296  0.4707 10.5597 
## 
## Number of Observations: 11156
## Number of Groups: 53
anova(m.hom.bas, m.het.mic.bas, m.het.subj.bas)
##                Model df    AIC    BIC logLik   Test L.Ratio p-value
## m.hom.bas          1  6 144451 144495 -72220                       
## m.het.mic.bas      2  7 144351 144402 -72168 1 vs 2   102.7  <.0001
## m.het.subj.bas     3 58 142054 142478 -70969 2 vs 3  2398.8  <.0001
## reactive

m.hom.rea <- lme(
  rt ~ trial.type * mic, 
  random = ~ 1 | subj,
  data = stroop.rt %>% filter(mic != "unknown"), 
  subset = session == "rea",
  control = cl1
  )

m.het.mic.rea <- update(m.hom.rea, . ~ ., weights = varIdent(form = ∼ 1 | mic))
m.het.subj.rea <- update(m.hom.rea, . ~ ., weights = varIdent(form = ∼ 1 | subj))

summary(m.hom.rea)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "rea" 
##      AIC    BIC logLik
##   165795 165839 -82891
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:       139.2    179.4
## 
## Fixed effects: rt ~ trial.type * mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          757.6     42.23 12469  17.937  0.0000
## trial.typei           79.9      7.29 12469  10.957  0.0000
## micmicro              54.0     47.33    52   1.140  0.2594
## trial.typei:micmicro  17.8      8.17 12469   2.175  0.0297
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.068              
## micmicro             -0.892  0.061       
## trial.typei:micmicro  0.061 -0.893 -0.068
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -4.6674 -0.5248 -0.1270  0.3272 11.4507 
## 
## Number of Observations: 12525
## Number of Groups: 54
summary(m.het.mic.rea)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "rea" 
##      AIC    BIC logLik
##   165469 165521 -82728
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:       139.3    167.6
## 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | mic 
##  Parameter estimates:
## micro fomri 
## 1.000 1.312 
## Fixed effects: rt ~ trial.type + mic + trial.type:mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          757.6     42.37 12469  17.881  0.0000
## trial.typei           79.9      8.94 12469   8.940  0.0000
## micmicro              54.0     47.44    52   1.137  0.2606
## trial.typei:micmicro  17.8      9.58 12469   1.856  0.0635
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.083              
## micmicro             -0.893  0.074       
## trial.typei:micmicro  0.078 -0.934 -0.080
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -4.9965 -0.5329 -0.1253  0.3398 12.0000 
## 
## Number of Observations: 12525
## Number of Groups: 54
summary(m.het.subj.rea)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "rea" 
##      AIC    BIC logLik
##   160559 160997 -80220
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:       138.2    174.8
## 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~1 | subj 
##  Parameter estimates:
##      115825      123117      130114      130518      132017      138837 
##      1.0000      2.5308      1.3131      0.4613      0.7545      1.1069 
##      141422      150423      158136      165032      171330      178243 
##      0.6145      0.8358      0.6950      0.7583      0.8991      1.4644 
##      178647      178950      197449      203418      204319      233326 
##      1.0869      0.9110      0.9354      0.6464      1.0887      0.8592 
##      300618      317332      346945      393550      448347      580650 
##      0.5252      1.0739      0.6571      0.7058      2.0377      1.7423 
##      594156      601127      765864      814649      843151      849971 
##      0.4918      1.1210      0.6392      1.0622      1.5554      1.1479 
##      873968      877168 DMCC1328342 DMCC1596165 DMCC1971064 DMCC2442951 
##      0.8275      1.5881      0.5940      1.7945      0.6483      0.3456 
## DMCC2803654 DMCC2834766 DMCC3062542 DMCC4191255 DMCC4854075 DMCC5009144 
##      0.6605      0.5432      0.8657      0.3361      0.3688      1.8097 
## DMCC5195268 DMCC5775387 DMCC6371570 DMCC6418065 DMCC6671683 DMCC6705371 
##      0.5287      1.0135      1.4468      0.6214      0.5124      0.4785 
## DMCC6721369 DMCC8033964 DMCC8050964 DMCC9441378 DMCC9478705 DMCC9953810 
##      0.4759      0.4892      0.7606      0.9560      0.8930      1.2308 
## Fixed effects: rt ~ trial.type + mic + trial.type:mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          760.5     41.95 12469  18.128  0.0000
## trial.typei           71.5      5.35 12469  13.364  0.0000
## micmicro              55.6     46.99    52   1.184  0.2417
## trial.typei:micmicro  12.5      5.86 12469   2.126  0.0336
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.050              
## micmicro             -0.893  0.045       
## trial.typei:micmicro  0.046 -0.913 -0.049
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -6.2200 -0.6222 -0.1659  0.4498 10.6857 
## 
## Number of Observations: 12525
## Number of Groups: 54
anova(m.hom.rea, m.het.mic.rea, m.het.subj.rea)
##                Model df    AIC    BIC logLik   Test L.Ratio p-value
## m.hom.rea          1  6 165795 165839 -82891                       
## m.het.mic.rea      2  7 165469 165521 -82728 1 vs 2     327  <.0001
## m.het.subj.rea     3 59 160559 160997 -80220 2 vs 3    5014  <.0001
plot(m.hom.pro)
plot(m.het.subj.pro)
m.het.subj.trial.type.pro <- update(m.hom.pro, . ~ ., weights = varIdent(form = ∼ trial.type | subj))

summary(m.het.subj.trial.type.pro)
## Linear mixed-effects model fit by REML
##  Data: stroop.rt %>% filter(mic != "unknown") 
##   Subset: session == "pro" 
##      AIC    BIC logLik
##   153344 153802 -76610
## 
## Random effects:
##  Formula: ~1 | subj
##         (Intercept) Residual
## StdDev:       121.8    159.9
## 
## Variance function:
##  Structure: Different standard deviations per stratum
##  Formula: ~trial.type | subj 
##  Parameter estimates:
##      115825      123117      130114      130518      132017      138837 
##      1.0000      1.3684      0.9483      0.7232      1.1539      1.0273 
##      141422      150423      158136      165032      171330      178243 
##      1.0303      1.3350      0.6193      1.6752      0.9424      1.5864 
##      178647      178950      197449      203418      204319      233326 
##      1.1050      0.7494      1.1987      0.7306      1.8857      2.0556 
##      300618      317332      346945      352738      393550      448347 
##      0.5359      0.6395      0.8655      1.6850      0.7852      1.1969 
##      580650      594156      601127      765864      814649      843151 
##      2.0928      0.6314      1.1921      0.5608      0.9475      1.1984 
##      849971      873968      877168 DMCC1328342 DMCC1596165 DMCC1971064 
##      1.8200      0.7821      1.3117      0.4071      1.6184      0.8052 
## DMCC2442951 DMCC2803654 DMCC2834766 DMCC3062542 DMCC4191255 DMCC4854075 
##      0.4382      0.6706      0.7042      0.9748      0.6480      0.6854 
## DMCC5009144 DMCC5195268 DMCC5775387 DMCC6371570 DMCC6418065 DMCC6661074 
##      0.5089      0.5855      1.3942      0.9953      0.5380      0.7065 
## DMCC6671683 DMCC6705371 DMCC6721369 DMCC7297690 DMCC8033964 DMCC8050964 
##      0.5219      0.7398      0.5582      0.7739      0.5176      0.8544 
## DMCC9441378 DMCC9478705 DMCC9953810 
##      1.4022      0.8540      0.9206 
## Fixed effects: rt ~ trial.type + mic + trial.type:mic 
##                      Value Std.Error    DF t-value p-value
## (Intercept)          735.5     35.44 11915  20.756  0.0000
## trial.typei           65.2      4.67 11915  13.940  0.0000
## micmicro              88.6     39.90    55   2.221  0.0305
## trial.typei:micmicro   8.7      5.45 11915   1.589  0.1121
##  Correlation: 
##                      (Intr) trl.ty micmcr
## trial.typei          -0.087              
## micmicro             -0.888  0.078       
## trial.typei:micmicro  0.075 -0.858 -0.091
## 
## Standardized Within-Group Residuals:
##     Min      Q1     Med      Q3     Max 
## -4.9122 -0.6275 -0.1734  0.4266 10.1056 
## 
## Number of Observations: 11974
## Number of Groups: 57
anova(m.het.subj.trial.type.pro, m.het.subj.pro)
##                           Model df    AIC    BIC logLik
## m.het.subj.trial.type.pro     1 62 153344 153802 -76610
## m.het.subj.pro                2 62 153344 153802 -76610